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This letter aims at resolving the issues raised in the recent short com- 
munication [I] and answered by 2 by proposing a systematic approxima- 
tion scheme based on non-mapped shape functions, which both allows to 
fully exploit the unique advantages of the smoothed finite element method 
(SFEM) [31 m E El El 111 [S] and resolve the existence, Unearity and posi- 
.^^ ^ tivity deficiencies pointed out in [l]. 

We show that Wachspress interpolants [TO] computed in the physical 
coordinate system are very well suited to the SFEM, especially when 
t/3 , elements are heavily distorted (obtuse interior angles). The proposed 

^^ ' approximation leads to results which are almost identical to those of the 

SFEM initially proposed in [3]. 

These results that the proposed approximation scheme forms a strong 
^ ' and rigorous basis for construction of smoothed finite element methods. 

I^l ' keywords: Smoothed Finite Element Method, boundary integration, Wach- 

r — , spress Interpolants, strain smoothing, rational basis finite elements, SFEM, 

CO ' isoparametric 

E ■ 1 INTRODUCTION 

The smoothed finite element method (SFEM) was first proposed in [5] . This new 
numerical method, based on gradient (strain) smoothing, is rooted in meshfree 
^1^ ■ stabilized conforming nodal integration |11| and was shown to provide a suite of 

Jh \ finite elements with a range of interesting properties. Those properties depend 

on the number of smoothing cells employed within each finite element (see [12j 
for a review of recent developments and properties) and include: 

• improved dual accuracy and superconvergence; 

• relative insensitivity to volumetric locking; 

• relative insensitivity to mesh distortion; 

• softer than the FEM. 

A rigorous theoretical framework was provided in [5l [8] and the method was 
extended to plates [9], shells [6] and coupled with the extended finite element 
method [H]. 



Table 1: Shape function value at different sites within an element (Figure ^) 
Site Node 1 Node 2 Node 3 Node 4 Description 
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0.0 


0.0 


Field node 
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0.0 


1.0 
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Field node 
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0.0 
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1.0 


0.0 


Field node 
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0.0 


0.0 


0.0 


1.0 


Field node 
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0.5 


0.5 


0.0 


0.0 


Side midpoint 
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0.0 


0.5 


0.5 


0.0 


Side midpoint 
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0.0 


0.0 


0.5 


0.5 


Side midpoint 
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0.5 


0.0 


0.0 


0.5 


Side midpoint 


9 


0.25 


0.25 


0.25 


0.25 


Intersection of two bimedians 





Figure 1: A four- node element divided into four smoothing cells. 

The essential feature of the SFEM is that no isoparametric mapping is re- 
quired, which implies that the approximation can be defined in the physical 
space directly, thereby providing freedom in the selection of the element geom- 
etry. 

In the initial paper [3] (Eq. (22) and reproduced here for simplicity Equation 
([TJ), non- mapped Lagrange shape functions are proposed as a possibility to 
calculate the shape functions at an arbitrary point within a smoothed finite 
element. It is then stated in the same paper (p863 last paragraph) that "unless 
state otherwise, we still use the averaged shape functions [3] for convenience." 
These shape functions are recalled in Tableland Figure ([T]), for ease of reading. 
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In fact, in our work on the SFEM [H |6l 111 HI [9] , and, to our knowledge, in 
all other work published to date fH HI [51 H] , these "averaged shape functions" 
have been used, with good results. 



Yet, [T] provides a critique of the SFEM stating that the approximation 
provided by Equation ^ are inadequate because: 

• they do not always exist (as described in the 1975 book [TU]): 

• they may not be positive everywhere in the element; 

• they may not be linear everywhere in the element. 

Because of this [T] disqualifies the current version of the SFEM and discredits 
the existing results of[3J|ll[51[SJ[71|51[n], despite the fact (also noted in [5]) that 
those non-mapped Lagrange shape functions of Equation ([1} were in general 
not used in the aforementioned papers. 

In this contribution, we show that it is possible to resolve the three issues 
mentioned by [I] about the Lagrange non-mapped shape functions while re- 
taining the advantageous features of the smoothed finite element method, in 
particular its ability to deal with extremely distorted meshes. 

2 WACHSPRESS INTERPOLANTS 

Wachspress [TU] presented a rigorous formulation for generating shape functions 
on arbitrary polygons, based on projective geometry. The Wachspress shape 
functions are unconventional compared to the polynomials used in the finite 
element literature, as they are in general rational functions, i.e, the ratio of two 
polynomials [13] . These shape functions have the following essential features of 
intcrpolants (for arbitrary n-sided polygons): 

• The Ni satisfy the partition of unity and Kronecker Delta properties; 

• Shape function, Ni is linear on sides adjacent to node i. 

• The Ni are linear complete. 

These properties make Wachspress intcrpolants effective to build approxima- 
tions on arbitrary n-gons, quadrilaterals in particular. Consider the quadrilat- 
eral (Q4), on which Wachspress interpolants will be formulated, shown in Figure 
^. Let hjhjh and I4 be the equations of the lines corresponding to each of 
the four sides of the quadrilateral, (1 — 2, 2 — 3, 3 — 4 & 4 — 1) respectively, and 
written in parametric form as: 

h (a;, y) = aix + bijj + ci ~ (2a) 

hix,!/) ^ a2X + b2y + C2 ^0 (2b) 

hix^y) ^ asx + bsy + cs^O (2c) 

hix, y) — 04a; -I- 642/ -I- C4 = (2d) 




Figure 2: A sample quadrilateral. 

where a.i,6j and ct for i =^ 1,2,3,4 are real constants. The wedge functions 
corresponding to each node, Wi, are defined so that they vary linearly along the 
edges adjacent to each node and vanish at the remaining nodes as [TO] : 



wi{x,y) = Ki hix.y) hix,y) 
W2{x,y) = K2 h{x,y) h{x,y) 
wz{x,y) = K3 h{x,y) li{x,y) 
Wi{x,y) = K4 li{x,y) l2{x,y) 



(3a) 
(3b) 
(3c) 
(3d) 



where the Ki are constants. In order that the Wachspress interpolants Ni, satisfy 
the partition of unity requirement, it is defined as: 
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(4) 



Remark 1 For a 4-sided polygon, the Wachspress rational basis interpolants 
degenerate to regular polynomials, i.e. '^Wi{x,y) is a constant. 

To illustrate the Wachspress interpolants for arbitrary quadrilaterals and to 
answer Zhang et al., [l]'s queries on 

• the existence of shape functions for arbitrary quadrilaterals; 

• positivity of the shape functions; 



linearity of the shape functions, 
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Figure 3: A parallelogram shaped element 



we choose the same parallelogram element as used in [Tj (of. page 1292, Figure 
2). This element is shown in Figure ^ for ease of reading. 
The shape functions for the parallelogram shaped element (see Figure ([3])) based 
on Wachspress interpolation write: 



Niix,y) = (y- 1) ( x- -y-1 



N2ix,y) = {l-y)[x- -y 



N3ix,y) ^y[x~-y 



N4x,y) = -y [ X - -y - I 



(5a) 

(5b) 
(5c) 
(5d) 



From the above equations, it can be easily verified that the shape functions 
derived using the Wachspress approach satisfy all the required properties (posi- 
tivity, Kronecker delta property, linear completeness) as shown in the book |10j . 
What is more, the quality of these functions is independent of the shape of the 
element, which make them ideal candidates for use in the context of the SFEM 
where highly distorted meshes are of interest given the absence of isoparametric 
mapping. 



For the mid-point of side 1-4, i.e., the point Q with coordinates (^j, 
Figure ([3])), the shape function values are: 
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as described in [I] (cf. page. 1293). 

Remark 2 By a systematic approach using rational basis functions, the shape 
functions on arbitrary n-gons, in particular Q4 's can be derived in the physical 
coordinate system, without the need for isoparametric mapping, and without the 
negative side effects of non-mapped Lagrange shape functions described in f^. 

Remark 3 The form of the Wachspress shape functions depend on the ele- 
ment geom,etry and hence have to be recomputed for each element in the mesh, 
which is computationally expensive. However, an important advantage is that 
their quality does not deteriorate for highly distorted elements, including con- 
cave domains. Interested readers are referred to \10[ \13\ \14^ for advances in this 
direction. 



3 NUMERICAL RESULTS 

3.1 Patch Test 

It is straightforward to check that the hnear patch test is satisfied down to 
machine precision. The interested reader is referred to the corresponding author 
to obtain a MATLAB code including this and other test cases. 

3.2 Bending of a cantilever beam 

As a second example, bending of a thick cantilever subjected to a parabolic load 
at the free end is examined as shown in Figure Q. The geometry is: length 
L=8, height D=4 and thickness t=l. The material properties are: Young's 
modulus E~3xl0^, and the parabolic shear force P=250. The exact solution 
of this problem is available in [15] . 

In this problem, two types of mesh are considered: one is uniform and regular 
and the other is an irregular mesh, with the coordinates of interior nodes given 
by: 



x' = X + (2rc — 1) air Ax 
y' = y+ (2rc - 1) air Ay 



(8a) 
(8b) 
(8c) 



Figure 4: Geometry and loads for a cantilever beam. 

where Tc is a random number between and 1.0, air G [0,0.5] is an irregularity 
factor that controls the shapes of the elements and Aa; and Ay are the initial 
regular element sizes in the x— and y— directions respectively. 

The domain is discretized with quadrilateral elements and two quadrilateral 
subcells for each element are used for the current study, denoted by SCkQ^j- 
Figure ([6]) shows a discretization of the cantilever beam with quadrilateral el- 
ements using an irregular mesh. A value of air = 0.5 is used to generate the 
irregular mesh. Under plane stress conditions and for a Poisson ratio i^~0.3, the 
exact strain energy is 0.0398333. Figure ([5]) illustrates the convergence of the 
strain energy and the rate of convergence in the energy norm of the elements 
built using the Wachspress interpolants compared with those of the SFEM re- 
sults (j8]) for regular quadrilateral meshes. The total strain energy and the 
rate of convergence for different irregular meshes (a^r = 0.2 and air = 0.5) are 
shown in Figures d?]), dH) and dH). 

It is seen that for regular meshes, the results are identical and for irregular 
meshes, there is a small difference, but the rate of convergence is not affected 
significantly. It is seen that with increase in mesh indejo, both methods converge 
to the exact solution. With the use of Wachspress interpolants, the question 
of positivity of the shape functions and existense of the shape functions for 
arbitrary quadrilaterals is resolved yet preserving the true essence of the SFEM. 

4 CONCLUSIONS 



This letter showed that it is possible to both retain the highly desirable features 
of the smoothed FEM (SFEM), and its true essence, i.e. strain smoothing and 
boundary integration, without sacrificing a rigorous approximation where the 
shape functions are known explicitly at any point of the smoothed finite element. 
The proposed interpolation scheme, known as Wachspress interpolation, pro- 
vides a suitable means to suppress the problems of definition, positivity and lin- 
earity associated with non-mapped Lagrange shape functions. It also provides a 
general framework to define approximation over arbitrary (possibly non-convex) 



^SCkQ4 implies k quadrilateral subcells for 4 noded quadrilateral elements 
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Figure 5: The convergence of the numerical energy to the exact energy and 
convergence rate in the energy norm with regular meshes for the cantilever 
beam problem: (a) strain energy and (b) the convergence rate. 
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Figure 6: Cantilever beam: (a) irregular mesh with extremely distorted elements 
and (b) zoomed view of the mesh. 
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Figure 7: The convergence of the numerical energy to the exact energy with 
irregular meshes {air = 0.2) for the cantilever beam problem. 

polytopes including curved edges or surfaces. 

The Wachspress approximation was conclusively tested within the context 
of the smoothed FEM, showing that it yields a method which passes the patch 
test and provides, even for very high element distortion, accurate and optimally 
convergent results in linear elasticity. 

The results also showed that utilizing the proposed Wachspress shape func- 
tions to build the approximation in the smoothed FEM leads to no significant 
difference in results compared to the widely used technique of "averaged shape 
functions" preconized in the original paper [3j and used in all the smoothed 
FEM work published to date [H El H il E [3 [6] . 

This finding reinforces the claim made in [7 that the approximation need 
not be known explicitly to build an optimal SFEM and dispels, also, the mis- 
conception of [1], whose detraction of the smoothed FEM rests on the misun- 
derstanding that this method requires the use of non-mapped Lagrange shape 
functions. 

It will be most interesting to investigate the behavior of the smoothed FEM, 
especially the choice of approximation, in the context of higher-order methods 
as well as non-polynomial enrichments |T21. 
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Figure 8: The convergence rate in the energy norm with irregular meshes (a^ 
0.2) for the cantilever beam problem: (a) SC2Q4 and (b) SC4Q4. 
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Figure 9: The convergence of the numerical energy to the exact energy and 
convergence rate in the energy norm with irregular meshes {air = 0.5) for the 
cantilever beam problem: (a) strain energy and (b) the convergence rate. 
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